Light propagation characteristics of turbulent plasma sheath surrounding the hypersonic aerocraft
1. IntroductionWith the development of hypersonic aircraft in near space, the plasma sheath surrounding the hypersonic aircraft has attracted the attention of many researchers. The plasma sheath is formed by the hypersonic movement of the aircraft.[1,2] When the aircraft moves at the speed of 10 Mach or more, the air around the aircraft begins to be compressed quickly and the temperature rises precipitously, generating a large amount of heat that ionizes the air. As the temperature rises, the ionization becomes intense and eventually forms a layer of plasma[3–5] which contains numerous free electrons, ions, and ablation particles. This layer is the so-called plasma sheath, which is inhomogeneous and variable with time. Meanwhile, recent studies show that the plasma sheath surrounding a hypersonic aircraft exhibits turbulent features.[6] For example, in 2014, Yang et al. studied the propagation characteristics of signals in time-varying plasmas. The results showed that the plasma parameters have time-varying characteristics due to turbulence effects of the plasma sheath or flight altitude adjustment and other factors.[7] In 2016, Shi et al. analyzed the influence of flight conditions (such as attack angle and small scale turbulence) on the transmission performance of the plasma medium. A new adaptive multistate Markov channel modeling method was proposed to present the dynamic effects of re-entry plasma sheaths on wireless channels.[8] In 2016, He et al. studied the dynamic effects of re-entry plasma sheath on electromagnetic wave propagation using the transmission matrix method. The established model included spatial turbulence and changes of electron density in plasma sheath. The amplitude and phase changes caused by the nonlinear characteristics of plasma were numerically simulated. It was found that the amplitude variation obeys logarithmic normal distribution, and the phase change obeys Rayleigh distribution.[9] Liu et al. found that the randomness of turbulent re-entry plasma sheaths can affect the propagation and scattering properties of electromagnetic waves. In terms of wave propagation, the randomness of turbulence caused amplitude and phase fluctuations, and the degree of fluctuation was closely related to altitude.[10] In 2017, Chen et al. used the finite difference time domain method to study the propagation characteristics of the terahertz wave in the plasma sheath of non-uniform space time. The study showed that the penetration of electromagnetic waves in the plasma sheath is closely related to temperature, pressure, and electron density.[11] In 2017, based on the hypersonic turbulence measurement of fractal theory, Li et al. derived the refractive index power spectrum index and the turbulent power spectrum. It was found that the turbulence in plasma sheath can seriously deteriorate the resolution of synthetic aperture radar imaging.[12] The existence of the turbulent plasma sheath has a great influence on the communication between the aircraft and radars on the land, and even results in communication blackout in some circumstances.[13–22]
To reveal the effect of plasma sheath on the communication of hypersonic aircraft, there have been many studies of the propagation of electromagnetic waves through the plasma sheath surrounding a hypersonic aircraft. In the early work, many researchers investigated the plasma sheath based on the aero-optical.[23–29] Within the aero-optical framework, the flow field was assumed to be low speed and the fluid was considered as air that is not ionized. In fact, air has been ionized into plasma in the hypersonic flow field and exhibits turbulent features. Therefore, researchers have attempted to solve hypersonic turbulence in aerospace engineering by employing other approaches. Li et al. studied the refractive index spectrum of plasma sheath and bit error rate of the laser through plasma sheath based on fractal theory.[30,31] Zhang et al. investigated the plasma channelʼs characteristics and performance.[32] In their work, the turbulent plasma sheath was treated as a fast fading wireless channel, and the simulation results showed that error floor took place for phase-shift keying and quadrature amplitude modulation, while frequency-shift keying with non-coherent detection was a promising way to mitigate the blackout problem. It is worth noting that they used low frequency waves. However, plasma sheath is a non-uniform medium with dispersion effect, and its dielectric properties are closely related to the frequency of the incident wave. High-frequency waves are more resistant to blackout effects than low frequency waves. Specifically, studies have shown that the transmission attenuation of terahertz wave is much smaller than that of millimeter wave in the same plasma.[33,34] However, to the best of our knowledge, the influence of electron density on refractive index has not been studied yet. Recently, we used the number density and Gladstone–Dale (G–D) constant of electrons to obtain the relationship between refractive index and total density of fluid instead of the G–D formula.[35] Further study showed that the electrons have a significant effect on the refractive index. Based on this conclusion, in the present work we investigate the propagation characteristics of high-frequency waves through turbulent plasma sheath surrounding the hypersonic aircraft.
The structure of this paper is as follows. In Section 2, the numerical calculation of the hypersonic flow field is given. Section 3 presents the relationship between refractive index and density in plasma sheath. Section 4 focuses on the influences of wavelength and electron on refractive index fluctuations. In Section 5, the ray tracing in the turbulent flow field is presented. Section 6 is a summary of the paper.
2. Modeling of plasma sheathIn the Cartesian coordinates, the three-dimensional Navier–Stokes (N–S) equations with the appropriate closure models, including thermochemical non-equilibrium, can be written as[36–40]
where
and
are conserved variable vector and source term vector, respectively,
and
are inviscid flux vector and the viscous flux vector, respectively (take the
x direction as an example), (
κ and
are thermal conductivity and turbulent thermal conductivity coefficient, respectively,
u,
v, and
w are velocity components along three directions,
k is turbulent kinetic energy,
is vibrational source term,
is the mass formation rate of individual species,
is the stress tensor component that can be closed by the turbulence model,
[41]
,
hi, and
are vibrational energy, sensible enthalpy, and density of each species, respectively, the subscript represents the species
,
ρ,
p, and
Ua are density, pressure, and average velocity of mixtures, respectively, and
Eint,
Et, and
Tint are internal energy, turbulence internal energy, and internal temperature, respectively.
Based on the finite volume method, the N–S equations can be solved by the advection upstream splitting method by pressure-based weight function[42] (AUSMPW+). Moreover, the implicit lower–upper symmetric Gauss–Seidel relaxation[43] (LU-SGS) is regard as time marching algorithm.
The commonly used turbulence two-equation models are k–(ε turbulence model, k–ω turbulence model, and shear stress transport (SST) turbulence model. To select a more suitable turbulence model, we use the electron density along the wall of the radio attenuation measurement (RAM) C-II flight experiment with a speed of 23.9 Mach and a height of 61 km. In the same flight scenario, based on the Navier–Stokes equations, different turbulence models are used to simulate the plasma sheath. The electron density along the wall of numerical simulation is compared with the RAM C-II flight experiment data, as shown in Fig. 1. The concept of goodness of fit in curve fitting is used to judge the simulation effects of different turbulence models, as shown in Table 1.
Table 1.
Table 1.
Table 1.
Simulation effects of different turbulence models.
.
Goodness turbulence model of fit test point |
k–ω model |
SST model |
k–(ε model |
First point |
0.1889 |
0.6156 |
0.6585 |
Second point |
0.9544 |
0.9650 |
0.9965 |
Third point |
0.9991 |
0.9997 |
0.9941 |
| Table 1.
Simulation effects of different turbulence models.
. |
From Table 1, we can see that the goodness of fit of different turbulence models is within the range of (0, 1) at each flight test point. The greater the goodness of fit, the better the simulation effect. Therefore, for the overall simulation effect, it is obvious that the simulation effect of k–(ε turbulence model is higher than that of the other two turbulence models. Therefore, k–(ε turbulence models are used to simulate turbulent plasma sheath.
It is worth noting that the k–(ε model, which is a two-equation model which employs partial differential equations to govern the transport of the turbulent energy, is suitable for high Reynolds number circumstances, where k and (ε are turbulence kinetic energy and dissipation rate, respectively.
The square root of k is regarded as velocity scale, and the length scale is represented as
The turbulence kinetic energy and dissipation rate can be written as
where (
μ,
, and
are dynamic viscous coefficient, turbulent viscous coefficient, and mass average velocity, respectively,
,
,
,
, and
P is defined as
,
ui,
uj, and
um are the time averaged velocity in tensor representation, and
is the Kronecker symbol.
In most cases, the air is partially ionized. So the gas mixture is usually considered as N2, O2, N, O, NO, NO+, and
. The reaction equation rate coefficients of the seven species are shown in Table 2. The forward rate coefficient kf and back rate coefficient kb are defined by the Arrhenius form of
,
.
Table 2.
Table 2.
| Table 2.
Chemical reaction equations and rate coefficients.[44]
. |
In Table 2, M represents collision, can be any one of N2, O2, N, O, and NO.
Next, we give the parameter of flow field as shown in Table 3. When choosing the turbulence model, we compare the numerical results with the real flight experiment, as shown in Fig. 1 and Table 1. At the same time, the correctness of simulation program is verified. The simulation program is then used to simulate turbulent plasma sheaths under different flight conditions as shown in Table 3. The Reynolds number Re is expressed as:
, where ρ is the density of fluid, v and L are the characteristic velocity and the characteristic length of the flow field, respectively, and (μ is the dynamic viscous coefficient. The Reynolds number can be used to judge whether the calculated state reaches turbulent state. Under hypersonic conditions, the Reynolds number is very large and the inertia has a greater influence on the flow field than viscous force. The fluid flow is also more unstable, and small changes of velocity are easy to develop and enhance, resulting in disordered and turbulent flow field. The calculated state in Table 3 satisfies the turbulent state. The exit adopts extrapolation interpolation conditions, and the wall surface is non-slip isothermal. The turbulent energy k is 0 in the wall surface. The dissipation rate is
, where
and
are the density and dynamic viscosity coefficient of flow, respectively.
Table 3.
Table 3.
Table 3.
Boundary values of different flight parameters.
.
Case No. |
U/m·s−1 |
Static pressure/N·m−2 |
Static temperature/K |
Turbulent kinetic energy k/m2·s−2 |
Dissipation rate (ε/J· kg−1·s−1 |
5017 |
5606.6 |
79.8 |
270.7 |
1571.7 |
3350014 |
5018 |
5936.4 |
79.8 |
270.7 |
1762.0 |
4210570 |
5019 |
6266.2 |
79.8 |
270.7 |
1963.2 |
5227154 |
4618 |
5895.4 |
131.3 |
266.9 |
1737.8 |
6911939 |
4818 |
5936.4 |
102.3 |
270.7 |
1762.0 |
5398940 |
| Table 3.
Boundary values of different flight parameters.
. |
3. Relationship between index of refractive and density of plasmaIn plasma sheaths, the mass of free electron is extremely small, and its contribution to mass density is negligible. However, it shows a special dominant effect for index of refraction in the ionized gas mixture.[39] The dependence of the index of refractive on electron density is completely different from that on neutral density. According to Alpher and Whiteʼs point of view,[45] the index-of-refractive difference
of electrons in the ionized gas mixture is at least 10 times greater than that of each atom. In addition, it is negative. Therefore, in the case of ionized gas, even if the degree of ionization is not large, its optical properties are governed by the free electrons. The total index of refractive of a partially ionized plasma is the sum of the index-of-refractive of electrons, positive ions, molecules, and atoms. According to plasma theory, the contribution of electrons to the index-of-refractive of the plasma can be determined as[46]
where
and
ω are the oscillation frequency of the plasma and the electromagnetic wave frequency, respectively. When
, equation (
9) can also be replaced by
[46]where the G–D constant of an electron gas is a negative value and is closely related to the wavelength.
is a function of the wavelength and the electron number density. Moreover,
ne is the phase refractive index of the electron gas. The positive sign is used when it is the group refractive index. When the incident light is monochromatic, the phase refractive index is required.
According to the Cauchy empirical formula, the refractive index of an atom is[47]
where both
A and
B are constants. Specifically,
that is, the refractive index of the atom is weakly correlated to the wavelength.
The contribution of positive ions to the refractive index is
where
C is a constant depending on what the substance is.
According to this analysis, the total refractive index difference of the plasma can be expressed as[47]
We can also get the G–D constant of each species in the mixture by means of the Cauchy empirical formula
where the values of
As and
Bs are shown in Table
4.
Table 4.
Table 4.
Table 4.
The As and Bs in Cauchy empirical formula.
.
Species |
As |
Bs |
Air |
2.227×10−4 |
5.67×10−15 |
N2 |
2.3251×10−4 |
7.70×10−15 |
O2 |
1.8654×10−4 |
5.07×10−15 |
NO |
2.1588×10−4 |
7.40×10−15 |
N |
3.2717×10−4 |
1.474×10−15 |
O |
1.8269×10−4 |
6.336×10−15 |
| Table 4.
The As and Bs in Cauchy empirical formula.
. |
4. Influence of wavelength and electron on refractive index fluctuationsBased on the numerical results given in Section 2 and the analysis of the refractive index of the plasma in Section 3, the refractive index fluctuation in turbulent plasma sheath is investigated in this section. In the following investigation, the light path is selected at the typical position of optical window of the aircraft. The special path is selected at the window of the aircraft and the three-dimensional coordinates of the path are x = [0.00674,0.0788], y = [−0.1555,−0.1225], and z = [0.0628,0.0527].
4.1. Effect of wavelength on refractive index fluctuationsFigure 2 shows the influence of electrons and the other six species on refractive index fluctuations at different wavelengths, by taking the flight parameters as 50 km and 18 Ma. As we can see, the influence of electrons on refractive index fluctuations is closely related to wavelength. The longer the wavelength, the larger the refractive index fluctuations. In contrast, the influence of the other six species on refractive index fluctuations is almost independent of wavelength.
4.2. Influence of different species on refractive index fluctuationsFigure 3 shows the influence of different species on the refractive index fluctuations at different wavelengths. We can easily find that the influence of electrons on refractive index fluctuations is different from that of the other six species. Moreover, the influence of electrons on the refractive index pulsation is closely related to the wavelength. The longer the wavelength, the more significant the influence of electrons on the refractive index fluctuation.
Figure 4 shows the distribution of electronic number density Ne at different flight parameters which are given in Table 5.
denotes the logarithm with the base of 10 for the electron number density. Comparing Fig. 4(a) with Fig. 4(b), we find that as the flight height increases, the electron density decreases. Comparing Fig. 4(b) with Fig. 4(c), we can see that as the flight speed increases, the electronic number density increases. Figure 5 shows the refractive index fluctuations along the propagation path at different flight parameters which are presented in Table 6. The results show that as the speed increases, the refractive index pulsation decreases. The larger the speed, the higher the degree of ionization. However, the G–D constant of electron is negative. Therefore, the larger the electron number density, the smaller the refractive index. As the height decreases, the refractive index pulsation decreases. In addition, the lower the height, the larger the density of air, and the more the electrons are ionized, which results in a smaller refractive index. In summary, refractive index pulsation, also called the “turbulent intensity” in the plasma sheath, is closely related to flight conditions and electron number density.
Table 5.
Table 5.
Table 5.
Flight parameters.
.
Case No. |
Altitude/km |
Mach No. |
Velocity/(m/s) |
4518 |
45 |
18 |
5864 |
5018 |
50 |
18 |
5936 |
5019 |
50 |
19 |
6266 |
| Table 5.
Flight parameters.
. |
Table 6.
Table 6.
Table 6.
Flight parameters.
.
Case No. |
Altitude/km |
Mach No. |
Velocity/(m/s) |
5017 |
50 |
17 |
5606 |
5018 |
50 |
18 |
5936 |
5019 |
50 |
19 |
6266 |
4618 |
46 |
18 |
5895 |
4818 |
48 |
18 |
5936 |
| Table 6.
Flight parameters.
. |
5. Ray tracing in a turbulent flow field5.1. Theory of ray tracingThe biggest difference between turbulent ray tracing and laminar ray tracing is that the refractive index interface must be taken in the corresponding delta neighborhood in turbulent ray tracing. According to the refractive index distribution of the flow field, the turbulent refractive index field of the blunt body is approximately circular refractive index field in the delta neighborhood. Based on the principle of the laminar ray tracing, a schematic diagram of the turbulent head and side ray tracing is established on the side of blunt body, as shown in Fig. 6.
The basic principle of turbulent ray tracing is similar to that of conventional laminar ray tracing. In the two-dimensional Cartesian coordinate system, the schematic of turbulent ray tracing is shown in Fig. 6. For the turbulent flow field, it is advisable to assume that the
th interface equation is
in the
th delta neighborhood. The tangent equation at the point
can then be written as
The angle between the tangent line and the
x axis is
In most cases,
, so the normal equation can be expressed as
Assuming that the expression of the
th incident ray is
and
, the incident angle at point
, which is the angle between the tangent line and the incident ray, is given by
According to the Lorenz–Lorentz formula and Snellʼs law, the relationship between the refraction angle and the incident angle can be given by the form
[47]
By using the inverse distance weighting algorithm, the refractive index n1, the densities
and
can be obtained. The interpolation formulas are as follows:
As shown in Fig.
6, the angle between the normal and
x-axis is
, so the
lith incident ray can be expressed as
Based on the method described above, the trajectory of the ray in the turbulent refractive index field can be obtained.
5.2. Simulation of ray tracingFigure 7 shows the changes of the angle between incident ray and the y-axis of the airborne coordinate along the propagation path. The flight parameters are chosen as 50 km and 18 Mach, and the initial incident angles are set to 35 degrees and 40 degrees. The wavelengths are
,
, and 1.55
m, respectively. From the figure, it is found that the angle between the incident ray and the y-axis of the airborne coordinate has two peaks. After further observation, we find that the first peak is located at the shock wave. Usually, there are shock waves in the hypersonic flow field. The density of flow field after the shock wave differs greatly from that before the shock wave, which results in the dramatic changes of refractive index changes and the first peak. In addition, since the density of flow field drops sharply away from the shock wave, the refractive index becomes larger. This is the reason why there is a second peak of the refraction angle. It can also be seen that the amplitude of the second peak is larger than that of the first peak. The longer the wavelength, the more obvious the difference.
Finally, figure 8 and 9 depict the changes of angle between the incident ray and the y-axis of the airborne coordinate along the propagation path with different wavelengths. The wavelengths of
and
are selected as incident waves. The fluctuation of the angle between the incident ray and the y-axis of the airborne coordinates at a wavelength of
is greater than that at
. As can be seen from Fig. 9, when the incident wavelength is
, there are approximately 10 degrees of offset at points A1 and A2 comparing with the initial incident angle in the case 5018. However, there are only about 0.12 degrees of offset with the wavelength of
at the points B1 and B2 in Fig. 8. This arises from the fact that the index of refraction is closely related to the wavelength. As discussed earlier, the longer the wavelength, the stronger the refractive index pulsation.
6. ConclusionIn this paper, based on the calculation and analysis of the turbulent hypersonic plasma flow field, the relationship between refractive index and density has been discussed. It is has been proven that among all the species, the electrons have a dominant effect on the refractive index fluctuations of the turbulent plasma sheath, and the refractive index fluctuation is also closely related to the wavelength and flight parameters. As flight speed increases or flight altitude decreases, the refractive index fluctuation becomes larger. In addition, the longer the incident wavelength, the larger the refractive index fluctuation. However, the G–D constant of electron is negative and that of the other six species is positive. Moreover, the influence of electrons on refractive index is much larger than that of the other six species. Therefore, the larger the density, the smaller the index of refraction of plasma sheath, which is exactly the opposite of aero-optical effect. Finally, the propagation of light in the turbulent plasma sheath is studied. When the incident wavelength is
, the maximum angle offset can reach about 10 degrees comparing with the initial incident angle in cases 5018 or 4818. The theoretical results obtained in this paper might be useful to study the communication of hypersonic aircraft.